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The resistance at the charge neutral (Dirac) point was shown by Checkelsky et al in Phys. Rev. 
B 79, 115434 (2009) to diverge upon the application of a strong magnetic field normal to graphene. 
We argue that this divergence is the signature for a Kekule instability of graphene, which is induced 
by the magnetic field. We show that the strong magnetic field does not remove the zero modes 
that bind a fraction of the electron around vortices in the Kekule dimerization pattern, and that 
quenched disorder present in the system makes it energetically possible to separate the fractional 
charges. These findings, altogether, indicate that graphene can sustain deconfined fractionalized 
electrons. 



I. INTRODUCTION 

The elementary excitations in the fractional quantum 
Hall (FQH) effect carry a fraction of the charge of the 
electron, as argued by Laughlin.EI Experimental evidence 
for the fractionally charged quasiparticles have been pre- 
sented for the V = l/i state, using shot noise measure- 
ments in Refs. ^ and ^ and by using a scanning tun- 
neling transistor in Ref. ^. It has remained an open 
question whether there are experimental systems in two- 
dimensions (2D) other than the FQH states for whi 
the elementary excitations carry a fractional chargel 
and for which the fundamental mechanism for electron 
fractionalization is different from that in the FQH effect. 

Recently, proposals for a mechanism to fractionalize 
the electron have been suggested for graphene-like sys- 
temsjj in which the electrons disperse according to the 
Dirac equation in 2D. The mechanism involves opening 
a mass gap in the Dirac equation, via the spontaneous 
breaking of a symmetry, the formation of a Kekule bond 
dimerization pattern. This mechanism is fundamentally 
different than the one in the FQH effect. Instead, it is 
closer concf!|Dtually to the mechanism for fractionaliza- 
tion in 1DJ3'eI which for 30 years was believed to be pecu- 
liar to ID and not possible in 2D. However, the strengths 
of the interactions in graphene are, alone, insufficient 
to lead to the Kekule instability. Moreover, even if the 
Kekule dimerization pattern formed, an axial gauge po- 
tential is needed to deconfine the fractionally charged 
quasmarticles that attach to vortices in the Kekule pat- 
tern.ta 

We are going to show in this paper that the presence 
of a magnetic field in graphene stabilizes a quasi-long- 
range order in the form of a Kekujp-distortion associated 
to a U(l) continuous symmetry,LHlil and opens an elec- 
tronic energy gap 2Ag. This quasi-long-range U(l) order 
can be destroyed by the unbinding of vortices due to ei- 
ther thermal- or disorder-induced fluctuations. We then 
show that fractional charges do attach to vortices in the 
Kekule pattern even in the presence of a strong magnetic 
field. Quenched disorder present in the system, remark- 
ably, makes it energetically possible to separate the frac- 



tional charges. These theoretical results on the formation 
of the Kekule pattern in the presence of a magnetic field 
and of static disorder explain the features observed, in. the 
experimental measurements by Checkelsky et a/.E3ii3 

Let us bricfiy discuss the known phenomenology of 
graphene in the presence of a magnetic field. Start- 
ing from the linearized Dirac spectrum £(fe) = ±hvp\k\ 
around the two non-equivalent Fermi points (valleys) 
single-particle Landau levels in a uni- 
form applied magnetic field B perpendicular to graphene 



(1.1) 
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integer, the cyclotron frequency 
/Ig, while the magnetic length is £^ = 
he/ {eB)x3 This single-particle spectrum leads to the 
quantization of the dimensionless Hall conductivity 

V = ha^y/e^ = 4 X (n + 1/2) = ±2, ±6, ±10, • • • , (1.2) 

if the Zeeman splitting and the Goulipmb or electron- 
phonon interactions are neglected JIjlJ as obsepjpd ex- 
perimentally for magnetic fields of up to 8 [T|£ili23 On 
the other hand, the observation in Refs. ^ and p4| of new 
plateaus at the filling fractions v ^ Q, ±1, ±4, for applied 
magnetic fields between 20 [T] and 40 [T] indicates that 
the 4- fold degeneracies of the Landau levels (1.1) must be 
lifted. For example, spin-splitting may be induced by the 
Zeeman or Coulomb interactions while valley-splitting 
{Kj^ and K_) may be induced by Coulomb or elcctron- 
phonon interactions. Many theoretical proposals to un- 
derstand these new plateaus and predict new ones have 
been made (for a review, see Ref. p5[ ). 

Here, we shall focus on the particularly interesting 
ri = Landau level. The dependence of the longitudinal 
resistance R^^ at this level (fixing the chemical potential 
/i = 0) for graphene deposited on a Si-Si02 substrate un- 
der intense uniform magnetic fields B of up to 32 [T] 
has been measured in Refs. 12 and 13, At fixed low 
temperatures (0.3-5[K]), R,j,^ grows exponentially with 
magnetic fields i? > 17 [T] upon approaching a sample- 
dependent critical field strength B^. The critical field 
B^ is larger for samples with lower zero- field mobility /i^. 
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There also appears to be a sealing regime with the seal- 
ing function R{b) = exp[2 x €3(1 — b)^^] fl/D with 
= 440, C2 = 1.54, b = B/B^{fj.J < 1, and g = 1/2 fit- 
ting well the magnetic field dependence of R^^ for sam- 
ples characterized by the mobility-dependent B^. This 
sealing is reminiscent of the Kosterlitz-Thouless (KT) di- 
vergence of the correlation length in the 2D classical XY 
model.o Finally, the temperature dependence of R^^ for 
a given sample with a fixed B very close to B^ approaches 
the thermally activated form exp(2Aj^/T) on the interval 
2K<T<16K.E3 

Our proposal to explain this KT scaling is closely re- 
lated to the scenario from Ref. However, the disorder 
was treated at the level of the self-consistent Born ap- 
proximation in Ref. ^ while the role of the Coulomb 
interaction, within the Hartree-Fock approximation, was 
emphasized. Here, we show that the eleetron-phonon in- 
teraction is critical to the selection of the Kekule insta- 
bility among all 15 instabilities that can open a gap at 
the Dirac point while preserving the electron charge as a 
good quantum number. lJ Moreover, we argue that a non- 
perturbative treatment of the disorder is also essential to 
explain the striking fact that B^[^^) is independent of T 
at low temperatures: it is the 2D classical random phase 
XY model that describes the phase transition observed in 
Ref. |l^. Within this classical random phase XY model, 
vortices unbind when the magnetic field is below _B^, and 
thus deconfined fractionally charged quasiparticles can 
be sustained in graphene. 

The line of arguments in this paper can be summarized 
as follows. For the relevant range of magnetic fields, the 
dominant energy scale originates from the kinetic energy 
modified by the orbital coupling to the magnetic field. 
The Coulomb interaction on the length scale of the mag- 
netic length - the long-range Coulomb interaction - is 
the leading subdominant energy scale. To these leading 
and subleading orders, the dynamics of interacting elec- 
trons in graphene preserves a combined U(4) symmetry 
arising fojuH the conservation of the electron charge, axial 
charge,E3Ej and two independent spin-1/2 and valley-1/2 
SU(2) symmetries. This continuous symmetry can only 
be broken spontaneously at zero temperature according 
to the Mermin- Wagner theorem. Any non- vanishing or- 
dering temperature thus requires the explicit breaking of 
this U(4) symmetry down to products of the U(l) sub- 
group enforcing the conservation of the electron quantum 
number with either finite subgroups or U(l) subgroups 
of SU(4). The bare Zeeman, eleetron-phonon, and short- 
range Coulomb interactions, although close to two orders 
of magnitude smaller than the long-range Coulomb in- 
teraction, provide these anisotropics. In particular, the 
eleetron-phonon interaction favors the Kekule instability 
over other instabilities, for example magnetic ones, so 
that, even if one could turn off electron-electron inter- 
actions, the Kekule instability would still compete with 
the Zeeman interaction. Hence, we first start by con- 
sidering the effects of lattice distortions in the ease of 
noninteraeting electrons and without the Zeeman inter- 



action in Sec. ||. In this approximation, we show that the 
Kekule instability triggered by the eleetron-phonon cou- 
pling can explain the experiments of Ref. [l3. In See. HI, 



we show that as long as the Kekule energy gap in the 
absence of the Zeeman coupling is greater than twice the 
Zeeman splitting (once that is turned on), the results of 
See. remain unchanged. However, for graphene, the 
value we obtain theoretically for the Kekule gap using 
the eleetron-phonon coupling alone is too close to this 
threshold of twice the Zeeman splitting (as reported in 
Ref. |2^). Therefore, a proper treatment of the electron- 
electron interactions is needed for selecting the correct 
low-temperature phase in Ref. ^ We restore the full 
Coulomb electron-electron interactions in Sec. IV. We 



first argue that the lattice nearest-neighbor Coulomb in- 
teraction boosts the eleetron-phonon coupling and sta- 
bilizes the Kekule order in spite of the competing Zee- 
man splitting. We then argue that the Kekul instability 
and the Zeeman splitting provide the needed anisotropics 
to reduce the symmetry and evade the Mermin- Wagner 
theorem that precludes ordering at a finite temperature. 
The finite ordering temperature is shown to depend on 
both the SU(4) symmetric electron-electron interactions 
and on the explicit anisotropics, including those from the 
eleetron-phonon interactions. We conclude that the re- 
sults obtained within the approximations of See. ^ re- 
main valid, although the temperature scale for the tran- 
sition is increased. 



II. 



THE KEKULE INSTABILITY 



In this section we are going to ignore the Zeeman and 
electron-electron interactions altogether. We will revisit 
this approximation in Sees. HI and [V, We shall first de- 
rive the mean-field Kekule instability induced by phonons 
in the presence of a magnetic field. We shall then prove 
the existence of zero energy states when the Kekule order 
supports a defect in the form of a vortex in the presence 
of a uniform magnetic field. 



A. The mean-field Kekule order. 

In the continuum limit, the orbital contribution from 
a magnetic field J5 = V A A to pristine graphene reads 



T-Lq = hvp ( k - — A ) • a 



(2.1a) 



where the 2D wave number k has the components fcj and 
and the Dirac matrices are 



•So. 



a, = (To 



(2.1b) 



The unit 2x2 matrices ctq; ''o; ^^'^ ^0 together with the 
Pauli matrices cr, r, and s act on the valley-1/2 {K^ 
and K_), sublattice-1/2 (A and B), and spin-1/2 (t and 
^) 2-dimensional subspaces of graphene, respectively. 
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The spectrum of Hq is 4-fold degenerate, for Hq com- 
mutes with the 16 Hermitean generators 



X, 



'To 
To 



^3 ' 



sj2.2a) 
^*J2.2b) 



(here k = 0,1,2,3) of the Lie group 

U(4) = U(l) X SU(4) ~ S0(6). (2.2c) 
Notice the electron charge U(l) subgroup generated by 



X, 



000: 



the spin-1/2 SU(2) subgroup generated by 



X, 



001' 



X, 



0021 



X, 



003' 



and the valley-1/2 SU(2) subgroup generated by 



X 



130' 



X. 



230' 



X 



300 



75- 



(2.3a) 



(2.3b) 



(2.3c) 



A Kekule instability is a periodic modulation of the 
nearest-neighbor hopping amplitude in graphene with the 
wave vec tor — K_ u In the continuum approxima- 
tion (2.1), it is represented by the Hermitean mass ma- 
trix 



(2.4) 



that] is parametrized by the real- valued numbers A^^ and 



lS.2>^ The Kekule mass matrices are 



'0' 



'0' 



(2.5) 



where a mass matrix is any Hermitean matrix 



X^ 



3 := c^Kj ® s„^, = 0,1,2,3, 

(2.6) 

that anticommutes with T-L^. [There are 16 mass matri- 
ces as shown in Ref. They generate a group U(4) 
distinct from the group U(4) whose generators are given 
by Eq. (2.2).] In the continuum approximation, only the 
8 generators 



OOk' 



3k. 



(2.7) 



(where k = 0, 1, 2, 3) of the Lie group ( ^.2| ) commute with 

•H:=-Ho+7Wk- (2.8) 

In the presence of the uniform magnetic field B = 
diA2 — d2Ai and the uniform Kekule order parameter 
A = A^+iA^ = Aoe^^o, for < Ag and < 6*0 < 27r the 
single-particle spectrum is the shifted Landau spectrum 



^0 



(2.9) 



where iV = 0, 1, 2, . • • . 

The mean-field value of the Kekule gap 2Aq induced 
by a magnetic field and by an electron-phonon coupling 



in the single-particle approximation (see also Ref. |^) IS 
obtained by balancing the gains in the electronic energy 
against the losses in the elastic energy of the lattice: 



(2.10a) 



where 



^i?elec=2x^ 



— e 



Nj 



B M- 



N=0 



2 X 



A 



is the electronic gain in energy, while 

2 /i 



2 



'A"0 



(2.10b) 



(2.10c) 



is the cost in elastic energy. Here, is the bond elastic 
constant and ties the gap 2Aq to the atomic displace- 
ment Uq for the Kekule order, which results because of 
the change in hopping matrix elements from the orbital 
overlaps. (The contribution to the electronic gain from 
the levels N = 1,2,- •• can be absorbed into a down- 
ward renormalization of the elastic rigidity K^/f^X^) 
We have assumed that the cyclotron energy is the largest 
energy scale of the problem. We have also introduced 
the area A and the number A/" = 4:A/ (3-\/3 a^) of sites of 



is t he lattice spacing). The abso- 



graphene (a ~ 1.42 [A 
lute minimum of Eq. ( 2.10 ) is obtained for the mean-field 
value 



^0 ~ MB 



fl A 



3%/3a^ /i 



(2.11) 



of the Kekule single-particle mass. If, following Ref. |3^, 
we make the estimate f\/K^ w 6.366 [eV], then the 
mean-field Kekule gap is estimated to be 



2 X Ao w 1.86 X B[K], 



(2.12) 



in units of Kelvin (with B measured in units of the Tesla). 
The single-particle electronic gap (2.12) overeptitaates by 
a factor of « 2 the measured activation gaps .1130 

We now discuss the all important phase fluctuations of 
the Kekule order parameter and the effective action that 
governs them. The complex- valued Kekule order param- 
eter A = Aq e'^ has a phase 9{r,t) that can fluctuate in 
space and in time. The dependence of the effective action 
on this phase stems from the dynamics of the phonons 
on the one-hand and from integrating out the electronic 
degrees of freedom on the other hand. 

In addition, we shall also account for the fluctuations 
in the hopping matrix elements due to disorder or to local 
curvature effects (tiny ripples) on the surface of graphene. 



6n 



75' 



(2.13) 
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The axial vector potential Ag encodes these disorder ef- 
fects in the hopping matrix elements. Notice the 



75 



(2.14) 



matrix, so that Ar^ couples with opposite signs to the two 
flavors at and K , in contrast to the electromagnetic 
vector potential A. _ 

As shown by Jackiw and Pijl3 this axial vector po- 
tential is intimately related to the phase of the Kekule 
order parameter. The axial gauge transformation 5* 

single-particle spinor for electrons, 
removes the phase of the Kekule order parameter, A — > 
Ae-'O ^ Ao, by shifting A^^^ ^ A^^ - 8^9/2 (where 
/i = 0, 1, 2). The effective action that follows upon inte-, 
grating the electrons can be inferred from symmetry.E2l 
The leading term is, at T = 0, [d^ = (vp^d^ , d^)] 



S 



J 

2h 



dtd'r{d^e~2A,^y 



(2.15) 



The stiffness J was computed at zero magnetic field in 
Ref. Now, we compute J in the presence of B. To this 
end, we consider the time independent Kekule texture 
with 9 = q - r and A^ = 0, or the gauge equivalent 9 = 



and Ar 



-9/2. 



We solve for the eigenvalues of H + 
5H treating A^ = ^ 9/2 to second order in perturbation 
theory (see details in Appendix iM . The electronic energy 
cost of = —q/2 to order jql determines the stiffness 



2tt 



(2.16) 



(This result coincides with the stiffness computed at zero 
field in Ref. once spin degeneracy is accounted there) . 

We shall treat the in Eq. (2.15) as a static random 
phase on phenomenological grounds dictated by symme- 
try. A microscopic justification for this step is the follow- 
ing. In a realistic sample, the matrix hopping elements 
between nearest-neighbor atoms will not be all the same, 
but there will be fluctuations caused by several factors, 
strain when in contact with the substrate being one fac- 
tor. The ran dom h opping amplitudes not only cause the 
phase in Eq. (2.15) to be random, but also the coupling J 
to be spatially varying. However, weak random exchange 
is irrelevant, and we thus simply focus on the uniform J 
case. 

The effective action on length scales larger than is 
thus that of the quantum random phase XY model in 
2D whereby we assume the white-noise distribution of 
vanishing mean (overlinc denotes disorder averaging) 

Aj^FjA^^gS.^S^'Hr-r'), i,j = l,2. (2.17) 

It is crucial to appreciate that g is a dimcnsionlcss cou- 
pling constant. It thus depends on the magnetic field 
through a function g{£/a,£g/£) of two arguments. The 
argument £/a is the mean free path associated with the 
disorder for the system at i? = in units of the lat- 
tice spacing a. The argument ig/£ must also be present 



because as decreases with increasing B, the single- 
particle overlaps due to disorder of (otherwise orthogo- 
nal) orbitals decrease. 

Next, we are going to argue tha t the 2 D ran dom phase 
XY model model defined by Eqs. (2.15)-(2.17) solves two 



experimental puzzles observed in Ref. ^ Experimen- 
tally, the resistivity of graphene is fit to a form consistent 
with the scaling of the correlation length ^ of the 2D clas- 
sical XY model, where the distance to the critical point 
is measured in terms of the reduced field b = B/B^. Why 
is it that the data fits a transition as a function of b in 
the 2D universality class of the classical as opposed to the 
quantum XY model? Why is it that the resistivities are 
rather temperature independent for temperatures differ- 
ing by a factor of 5, from 0.3 K to 1.5 K? 

First, the measured KT scaling could be explained if 
the fluctuations driving the transition were ID but of 
quantum origin. For instance, if the transition had to 
do with edge physics. This idea has been proposed in 
Ref. The problem with this scenario is that it re- 
quires 1) a large Zceman coupling, 2) a critical density 
of magnetic impurities, and 3) the ordered phase in a ID 
quantum system does not survive finite temperatures so 
that the size of the region in which the correlation length 
scales in a KT fashion near the quantum critical point is 
bounded from above by an exponentially small crossover 
temperature. 

Second, a quantum critical point of the 2D quantum 
XY model ( 2.15| ) without static random phases does not 
obey the KT scaling that fits the data in Ref. |l3| over 
3 decades of values of longitudinal resistances R{b). On 
the other hand, a KT scaling is compatible with the clean 
limit of model (2.15) if the relevant range of temperatures 



for which quasi-long-range order holds is large enough so 
that quantum critical fluctuations can be neglected (see 
Ref. p3| ). However, even if the quantum critical fluctu- 
ations can be neglected when the temperature changes 
from 0.3 K to 1.5 K, the position of the critical field B^ 
should move with temperature, for the location of the 
clean KT transition is where J{B) = 2T/tt. 

Third, the disorder in Eq. ( pl5| ) can explain the ob- 
served KT scaling. Indeed, the 2D classical random phase 
XY model has the remarkable feature that its boundary 
separating the quasi-long-range ordered phase from the 
paramagnetic phase runs (nearly) parallel to the temper- 
ature axis when T < rj^rp/2. where ^kt Tf 7/2 (This 
feature is non-perturbative in disorder .l3c3) However, no 
consensus has emerged on the exponent g that character- 
izes the diverging correlation length exp(6/a;^) upon ap- 
proaching the boundary (it is debated if it is g = 1/2 or 
1). The phase diagram is shown in Fig. |l|. Also shown are 
the trajectories of coordinates {T/ J{B), g{£/a,ig/£) ) 
as the magnetic field B is changed for fixed £ and 
T. When these trajectories enter the quasi-long-ordcrcd 
phase through the disorder-controlled boundary parallel 
to the T/ J axis, one has the condition 

g{£/a, £b/£) = tt/S, i.e., B = B^{a, £/a) (2.18a) 
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Increasing B for a fixed T 




FIG. 1: (Color online) Phase diagram of the 2D classical 
random phase XY model after Ref. |3^. The horizontal axis 
represents the temperature T in units of the stiffness J. The 
vertical axis represents the dimensionless disorder strength 
g (the variance for white-noise and Gaussian correlated ran- 
dom phases). The shaded region depicts the quasi- long-range 
ordered phase for which the interaction between vortices is 
logarithmic. The complementary white region depicts the 
paramagnetic phase for which the interaction between vor- 
tices is screened beyond the screening length ^. At the tran- 
sition line that separates the two phases, vortices undergo a 
confining-deconfining transition. We interpret the divergence 
of the resistance of graphene subjected to a strong increasing 
magnetic field B measured in Ref. |l^ as being governed by 
the increase of the screening length (,{B) in the paramagnetic 
region along a RG trajectory parametrized by B. Here, T/J 
and g are both decreasing functions of B. 



which is temperature independent as long as 
T < rKT/2, i.e., T < Ao(i?J/8. 



(2.18b) 



One can estimate the upper range of temperatures for 
which the transition is KT-like as a function of the scahng 
parameter b = B/B^ use d in Ref. fsl For B, sa 29.1 



[T], one finds using Eqs. ( |2.12D and ( p^lSbl) that the 2D 
classical scaling as a function of b holds as long as T < 3.5 
K. 

These results thus explain rather well the features ob- 
served in Ref. |l^. Accordingly, the Kekule distortion, 
which does not appear to occur in graphene at zero mag- 
netic field, is induced by the presence of the magnetic 
field above a critical field B^ that depends on the disor- 
der strength and is thus sample dependent. 



B. Zero modes in a uniform magnetic field 

One of the surprising features in graphene is that vor- 
tices in the Kekule dimerizatipn pattern bind fractional- 
ized excitations when B = OU We now show analytically 
that fractionalization is robust to the presence of B, i.e., 
that there is a zero mode due to vortices that binds a 
fraction of the electron for any value of B. 



We thus seek e = solutions of the Dirac operator ( |2.8| ) 
when the Kekule order parameter supports the vortex 



(2.19) 



To this end, we use polar coordinates and the symmetric 
gauge A = B{—y,x)/2. Without loss of generality, we 
restrict the vorticities to tt-v = ±1. It follows from the 
sublattice symmetry in the problem that, as in the case 
when B = 0, the zero- mode solutions have support on 
either the A (for riv = — 1) or B sublattice (for — +1). 
The solutions in thepresence of the magnetic field (see 
details in Appendix p[) are 



*A(r) 





where 



u{p) 
v{p) 



e+^<'o/^D_^{p/e^) , 



(2.20) 



(2.21a) 
(2.21b) 



Dp{z) is the parabolic cylinder function,^ and the di- 
mensionless ratio 



2A§ 



(2.21c) 



The existence of a single midgap state for a vortex 
(riv = +1) or anti- vortex (n^ = — 1) implies that thesg 
topological defects bind half of the electron charge ±1/2.0 
We have verified numerically (see Fig. ||) that the effect of 
breaking the sublattice symmetry with a staggered chem- 
ical potential p^ (a coexisting charge density wave) is to 
move the midgap states away from p ~ 0. According to 
Fig. ^ the fractional charge bound to a vortex acquires 
a dependence on the value of p^ . Remarkably, we can fit 
this dependence in the presence of the strong magnetic 
field B with the results obtained in Refs. ^ pSl |39|, 



and 29 in the absence of the magnetic field! Altogether, 



our analytical and numerical results show that the frac- 
tionalization of the electric charge of electrons is robust 
to any magnetic field. 



III. 



COMPETING ZEEMAN ENERGY 



So far, we have ignored the Zeeman energy. The reason 
to do so is two-fold. First, the Kekule order parameter is 
a spin singlet and thus all preceding considerations apply 
in the presence of a Zeeman interaction, if it can be estab- 
lished that the Kekule instability survives the presence 
of the Zeeman interaction. Second, as we show below. 



the mean-field Kekule instability ( 2.11 ) is unchanged by 
the addition to H defined in Eq. (2.J) of the Zeeman 
interaction 



A2 CTq ' 



(3.1) 
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Q(e) 




FIG. 2: Fractionalized charge Q in units of the electric charge 
e of electrons in the presence of a magnetic field B = 20 
(T) as a function of the ratio between the staggered chemical 
potential ^is{oo) and the Kekule gap Ao(c)o) very far away 
from a single vortex (filled circles) or very far away from a 
pair of topological defects made of a vortex and a half axial 
flux (open circles). The lattice model consists of 113 x 80 
sites. The solid lines are the analytical results presented in 
Ref. |2S| when B = 0. 



if 

while there is no Kekule instability if 

An 



< —0 < A, 



the case of 



(3.2) 



(3.3) 



(3.4) 



being a me an-field first-order transition. Hence, if con- 
dition (3^) holds, it is legitimate to ignore the Zeeman 
energy altogether when computing the magnitude of the 
Kekule gap. 

This is not to say that the Zeeman interaction does not 
play an essential role, for it is it together with the Kekule 
order parameter that establishes the desired pattern of 
explicit symmetry breaking 

Ho ^0 +'Wk +^z, (3-5a) 
U(4) ^ U(l) X U(l) X U(l) X U(l) (3.5b) 

with the residual unbroken U(l) symmetry generators 



"'^000' ^003' A]_Xj_3o A2X23Q, Aj_Xj 
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A2-''^233' 



(3.5c) 

respectively. 

To understand results (3.2)-(3.4), it suffices to observe 
that the presence of the Zeeman interaction (3.1) changes 



the single-particle spectrum (|2.9|) to 



'-N.± 









+ A2, 






^A2 


-Az, 






^A2 


+ A2, 



(3.6) 



where N = 0,1, 2, ■ ■ ■ . I t then follows that the gain in 
the kinetic energy ( 2.10b ) is changed to 



SE. 



clcc 



2 X 



4 X 



A 



2nel 
A 



e(Ao-Az)(Ao-Az) 



^ N=l 



E 



(3.7) 



-N) 



where Q is the Heaviside function. The difference be- 
tween the total (electronic and phonon) energy with and 
without the Kekule order is 



SE 



-2 X 



A 

2nel 



e(A„-A2)(A„-A2) + ^^A^ 



(3.8) 

according to Eq. (2.10c). Condition (3.2) follows from 
seeking the global minimum of SE{Aq) with Aq > and 
5E{A,)<0. ^ . rn . 

The measurements from Jiang et al. in Rcf . ^ give the 
estimate 



2 X Az w 1.3 X S [K] 



(3.9) 



for the Zeeman splitting. Hence, although the phonon- 



induced Kekule gap (2.12) is larger that the o ne m easured 
in Rcf. E4, it is short of satisfying condition (3.2), unless 



it is enhanced by another mechanism, say by the short- 
range part of the Coulomb electron-electron repulsions. 
Indeed, there are contributions other than phonons that 
also favor the Kekule instability. For example, in Rcf. |^, 
the nearest-neighbor (extended Hubbard) potential was 
found to favor the Kekule instability for strong enough 
coupling. In the next section, we turn to the discussion of 
the Coulomb interaction, and the issues of its symmetries 
and anisotropics. We will show that the anisotropics are 
essential to lower the symmetry of the problem to either a 
discrete symmetry or an Abelian U(l) symmetry, so that 
ordering can take place in the 2D system. The fact that a 
Kosterlitz-Thouless transition is seen in the experiments 
of Ref. |lj suggests that the interactions arc sufficient to 
boost the U(l) Kekule ordering over the Zeeman term. 



IV. ELECTRON-ELECTRON INTERACTIONS 

We have ignored Coulomb electron-electron repulsions 
so far. In the absence of magnetic fields and by power 
counting, short-range interactions are marginally irrel- 
evant while the 3D long-range Coulomb interaction is 
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marginal. In the latter case, it is an empirical fact that 
the Coulomb long-range interaction is marginally irrele- 
vant in the samples of graphene studied to this date. 
The 3D long-range Coulomb interaction takes the form 



1 



Cb 



- d^r dV'p(r)Fcb(r)/5(r'). (4.1) 



The local electronic density operator 
p(r) := (1p^) (r) 



(4.2) 



is constructed from the 8-component spinor-valued 
operators T^'^{r) and ^/'(r) that follow upon second- 
quantization of the Dirac single-particle 8x8 Hamilto- 
nian (2.1): 



H, 







dVV^(r-)-Ho(r)V'(r). 



(4.3) 



Finally, the 3D Coulomb two-body repulsive interaction 
is 



e 



(4.4) 



with e the effective dielectric constant resulting from the 
substrate below and the air above graph ene. 

I shares the 



The long-range Coulomb interaction (4.1 



global U(8) symmetry of the coarse-grain electronic den- 
sity (4.2) under 



(4.5) 



This global U(8) symmetry is broken down t o its U(4) 
s ubgr oup generated by the generators in Eqs. (2.2a) and 
( 2.2b ) that leave the kinetic energy in Eq.(O) invariant. 

The characteristic energy scale of the long-range 
Coulomb interaction in the presence of a magnetic field 
is estimated in Refs. |l|, and ^ to be 



120/BfKl 



e 



(4.6) 



in units of Kelvin (with B measured in units of the Tesla). 
This makes the long-range Coulomb interaction the sec- 
ond largest energy scale for the range of magnetic fields 
relevant to Ref. O after that arising from the kinetic en- 
ergy (p|). 



The U(8) symmetry (4.5) of the long-range Coulomb 
interaction is not exact. On the length scale of the hon- 
eycomb lattice spacing, the local electronic density oper- 
ator is not given by the coarse- grain approximation (4.2) 
but it resolves the sublattice-1/2 and valley-1/2 quantum 
numbers down to the point group symmetry of the tight- 
binding model. Consequently, the on-site or the nearest- 
neighbor density-density repulsive interactions break the 
U(8) symmetry down to the product between the U(2) 
subgroup that generates the conservation of the charge 
and of the spin-1/2 quantum numbers on the one hand 



and of the discrete point group on the other hand. The 
spin-1/2 SU(2) symmetry is in turn broken by the Zee- 
man interaction ( |3.lD . 

The characteristic energy scale for the Coulomb repul- 
sion on the scale of the lattice spacing must be of the 
order 



X fia/iB) 



(4.7) 



with / a dimensionless function obeying /(x — > 0) = 0. 
If the asymptotic expansion /(x — 0) ~ a; holds, it then 
follows that the characteristic energy scale of the lattice 
Coulomb interaction grows linearly with B 



J- X fia/is) 



B 



(4.8) 



and could thus be comparable to the Zeeman split- 



ting (3.9) or to the phonon-induced Kekulc gap (2.12). 
We shall assume this to be the case. 

Because of the Mermin- Wagner theorem, the long- 
range Coulomb interaction cannot break spontaneously 
the combined U(4) spin and valley symmetry of the ki- 



netic energy "Hg from Eq. (4.3) at any non- vanishing 
temperature. However, any subdominant interaction 
that breaks the U(4) symmetry down to a U(l) or a fi- 
nite subgroup of U(4) (say its center Z4) can establish 
a non-vanishing-temperature transition to a quasi-long- 
range or long-range-order phase, respectively. Hence, 
any of the Zeeman, the phonon-electron, or the short- 
range Coulomb interactions can induce a non- vanishing- 
temperature phase transition. This brings about two 
questions. What is the leading instability and at what 
temperature does the leading instability takes place? 

First, since any finite-temperature instability must be 
associated to a discrete symmetry or to a U(l) sym- 
metry subgroup of U(4), the empirical observation of a 
KT-like transition rules out any discrete symmetry sub- 
group. This empirical fact can be understood at the 
mean-field level from three observations, (i) We showed 
in Sees. ^ and [II that phonons favor a Kekule instabil- 
ity with a gap larger (but not substantially larger) than 
the Zeeman splitting, (ii) Moreover, a repulsive nearest- 
neighbor Coulombj-jinteraction reinforces the phonon- 
induced instability.Ll (iii) An exhaustive list of instabil- 
ities in graphene can be found in Ref. Because the 
cyclotron energy is the dominant energy scale for the rel- 
evant range of magnetic field, we can ignore U(l) insta- 
bilities associated to superconductivity. This still leaves 
several competing U(l) order parameters (see Table 1 in 
Ref. |2^) but of these only one, the Kekule instability, is 

also favored by phonons. 

Second, we have seen in Eq. ( |3.6| ) that the Zeeman in- 
teraction opens a single-particle gap at half-filling in 
the Landau level v = 0, i.e., the Slater determinant made 
of al l ly = single-particle energy eigenstates of Hq in 
Eq. ( ^.l[ ) with negative eigenenergies has a larger energy 
than that of Hq + T-L^, ^1 half-filling. The Zeeman term 
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is also a symmetry bre akin g field that lowers the U(4) 
symmetry of Hq in Eq. (2.1). A ny o f the 15 traceless gen- 
erators of the U(4) symmetry (2.2) would do the same, 
or, more generally any 15-component uniform real-valued 
"SU(4)-symmetry-breaking magnetic field" b = (b^) in 



15 



n{b) ■.= Y.b^T^ 



(4.9a) 



a=l 

^ Xooi,... ,Ti5:= 4^X303, (4.9b) 



16 



16 



opens a single-particle gap if added to Hq. [The factor of 



l/vl6 insures the normalization tr T^^T,^ = for the 

15 generators of SU(4).] 

These single-particle considerations are also valid in 
the presence of the Coulomb interaction. As is the case in. 
double-layer non-relativistic quantum Hall systems,EJC3 
the long-range Coulomb interaction in the v — Q Lan- 
dau level is minimized by breaking spontaneously the 
SU(4) symmetry in Eq. (2.2) in the jS^ipf way as a quan- 
tum SU(4) ferromagnet would do,c2rcj i.e., in such a 
way that ip'' {r)7i{b)ip{r) acquires spontaneously a finite 
expectation value independent of r in the many-body 
ground state for some "SU(4)-symmetry-breaking mag- 
netic field" b = (b^). Such a many-body ground state 
is symmetric under any exchange of a pair of its SU(4) 
indices. Consequently, the Pauli principle can maximize 
the spatial separation of electrons to lower the Coulomb 
interaction without penalizing the kinetic energy due to 
the magnetic field. 

If the Coulomb interaction is treated in a mean-field 
approximation by which Hi^^ — > 7i{b) for some "SU(4)- 
symmetry-breaking magnetic field" b ~ (^a)i ^he pres- 
ence of a single-particle mean-field gap justifies integrat- 
ing the fermions in a gradient expansion. If this is done at 
any small but non- vanishing temperature, we can ignore 
quantum fluctuations to a flrst approximation. There 
follows an effective classical theory that describes how 
small thermal fluctuations about the SU(4)-symmetry- 
broken ground state restore the full SU(4) symmetry, as 
dictated by the Mcrmin- Wagner theorem (we always as- 
sume that the dominant cyclotron energy prevents any 
superconducting instability). On symmetry grounds, this 
classical theory is a non-linear-sigma model (NLSM) on 
the target manifold G/H given by the N — 1-dimensional 
complex projective space 



CP 



U(iV)/U(iV - l)xU(l) 
SU(iV)/S(U(iV- l)xU(l)) 



(4.10) 



with iV = 4 (see Rcfs. |6|, EtI, and || for the introduc- 
tion and overview of CP^ NLSM). The little group H 
is here the largest subgroup of G that leaves the ground 
state ip{b) of H{b) invariant. It is made of the direct 
product of the subgroup SU(3) spanned by all unitary 
transformation in the subspace orthogonal to the direc- 
tion of Tpib) in the combined spin and valley subspace and 
the subgroup U(l) spanned by all rotations about tp{b) in 



the very same spin and valley subspace. If we also include 
the symmetry-breaking subleading interactions (the Zee- 
man, electron-phonon, and electron-electron short-range 
Coulomb interactions), this effective theory depends on 
two dimensionless couplings. There is the reduced tem- 
perature t. There also is the coupling A <C t that breaks 
the SU(4) symmetry to subleading order so as to enforce 
the desired pattern (3.5) of explicit symmetry breaking. 
By power counting, t is marginal while A is strongly rel- 
evant. (These couplings have already been renormalized 
by quantum fluctuations.) 

On the one hand, if we ignore the anisotropy A, the re- 
duced temperature t flows away from the unstable infra- 
read fixed point t = to strong coupling with the beta 
function (see Refs. pQ pU p3, and 53) 



dt 



(4.11) 



n=l,2,3, 



that vanishes to order n = 1 and is given by the quadratic 
Casimir in the adjoint representation C^, 



(32 = kx C,„ = kx N 



(4.12) 



with iV = 4 to one loop order (n = 2), up to the multi- 
plicative factor k that depends on the convention made 
for the normalization of the reduced temperature (k = 1 
in Rcfs. 11^, |5^, |5^, and ^3|). Correspondingly, as t grows 
under this renormalization flow the correlation length 
^{t) decreases. 



^(t) ~ oexp 



l/(fc X N) 



t 



(4.13) 



with = 4. On the other hand, the mass scale \^^'^ 
grows away from the Gaussian (unstable) fixed point t = 
A = to leading order in the flow of the renormalization 
group according to 



/3a(A) 



9^ n^ 

a— = 2A. 
oa 



(4.14) 



The symmetry crossojiGrjs estimated from the condition 
(^/a)-2 « A< 1, i.e.^ita 



2/ikxN) 
ln(c/A) ^ 



(4.15) 



with A ^ c a positive constant of order one fixed by the 
short-distance physics. Once Eq. ( 4.15 ) is fulfilled, it is 
the classical renormalization group flow in the unbroken 
subgroup of U(l) of SU(4) that takes over. 

We conclude that the Kekule mechanism can explain 
the experiment of Ref. ^ if (i) the phonon-induced 
Kekule instability can be boosted by the nearest-neighbor 
short-range Coulomb interacti on t hat favors the Kekule 
instability so as to satisfy Eq. (3.2), (ii) the renormalized 
SU(4)-symmetric energy scale >/gu(4) [^^ opposed to the 
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bare estimate (4.6)] that enters the reduced temperature 
in the NLSM satisfies 



l/(fc X 4) 



X 0.12 X Jgu(4) [K] 



(4.16) 



since A is of order a/£g ^ 1/50 for the r elevant range of 
magnetic fields and if c = 1 in Eq. (4.15). This estimate 
will be lowered by quantum fluctuations and by disorder 
(as described in Sec. |l|for the disorder). 



V. DISCUSSION 

In summary, we have argued that sufficiently large 
magnetic fields and sufficiently small temperatures (in- 
side the shaded region of the phase diagram of Fig. |^) 
stabilize the Kekule quasi-long-range order in graphene. 
This quasi-long-range order can be destroyed by the un- 
binding of vortices due to either thermal- or disorder- 
induced fluctuations. The role of the magnetic field B 
is here 2-fold. First, B enhances the density of states at 
the Dirac point, thereby favoring the Kekule distortion 
with the help of phonons. Second, a large B prevents the 
disorder from filling a single-particle gap by suppressing 
the Landau degeneracy. The transition belongs to the 2D 
classical random phase XY model and is characterized 
by the deconfinement vortices for sufficiently small mag- 
netic fields or sufficiently large temperatures (outside the 
shaded region of the phase diagram of Fig. I). Each of 
these vortices, as we have shown here, bind a fraction of 
the electron charge for any magnetic field strength. 
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Appendix A: Stiffness 

In this appendix, we start from Hamiltonian ( |2.8| ) with 
the uniform magnetic field B ~ dAy/d^ — dA^/dy and 
with the uniform Kekule order parameter A = AQg'^o. 
We choose to work in the Landau gauge A = {^By, 0). 
We are after the Kekule phase stiffness. 

To this end, we first need the eigenfunctions and 
eigenenergies of Hamilto nian (2.8). That the eigencn- 
ergies are given by Eq. ( |2.9| ) follows from the fact that 
the uniform Kekule order parameter anticommutes with 



the covariant derivative. Hence, any positive eigenen- 
ergy of is nothing but the sum of two positive terms 
adding in quadrature. One is the eigenvalue of the square 
of the covariant derivatives in 



Eq. (^ 



The other is 

the eigenvalue of the s quar e of the uniform Kekule or- 
der parameter in Eq. (2.8). To obtain the (degener- 



ate) eigenfunctions of any Landau level (2.9), we antic- 
ipate that any eigcnfunction of Eq. (2.8) factorizes into 
the x-dependent phase e**^^ times a 4-component spinor 
that depends solely on y through the dimensionless co- 
ordinate ^ = {y /£g)+k£g. As it does for the conventional 
(non-relativistics) Landau levels, the wave number k en- 



t) that scales with 
frelativistic') source 



codes a degeneracy of the spectrum (2.' 
the area of the system. An additiona 
of degeneracy of the spectrum (2.9) arises because of the 
dimensionality 4 of the Dirac matrices. This finite de- 
generacy in the sole presence of a uniform magnetic field 
is selectively lifted by the uniform Kekule order param- 
eter. Indeed, the single-particle eigenstates with TV = 
are non-degenerate for each wave number fc as a result of 
the Kekule instability and given by 



fe,0,± 



„ikx 



V2 








(Al) 



On the other hand, for any wave number k, positive in- 
teger iV = 1, 2, • • • , and sign ±, the pair of ortho normal 
eigenfunctions 



k,N,± 



„ikx 



V2 



and 



k,N 











(A2) 



(A3) 



remain 2-fold degenerate in spite of the Kekule instabil- 
ity. Here, (pj^ with = 0, 1, 2, . . . are the orthonormal 
eigenfunctions of the one-dimensional quantum harmonic 
oscillator. 

Next, the Kekule phase stiffness can be computed by 
adding all t he c hanges in the negative single-particle en- 
ergy levels (^^) up to the Fermi level = that are 
induced by the twist A — )• Acxp(i<7 • r) in the Kekule 
order parameter. 

Before doing this, it is convenient to gauge out the 
spatial dependence of the Kekule order parameter with 
the help of the pure axial gauge transformation = 1/"^, 
where U = e""' '^'''^/^. Under this transformation, the 
Hamiltonian becomes "H' = T-l + V, where Ji is given 
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by ( p.8[ ) with the Landau spectrum ( p.9[ ) and V = Up q • 
a 75/2. 

We treat V as a perturbation, assume that /iWp|<7| 
is smah compared to Aq, and compute the change in 
the spectrum of % induced by V up to second order in 
degenerate perturbation theory. Because the perturba- 
tion V does not couple eigenstates with different wave 
number k, the degeneracy of each unperturbed Landau 
level ( |2.9D is conserved. Furthermore, these second-order 
shifts e)^\^ _|_ cancel pairwise for all 2-fold degenerate lev- 



els: 



^i'^TV + + Sfc^Tv _ = when N = 1, 2, • • • . Hence at 
half- filling (/^ = 0), the second order shifts from all states 
with N = 1, 2, • • • cancel pairwise, leaving the unpaired 
iV = states as sole contributors to the total energy 
shift: 



5E = 2x^e 



(2) 
k,0, 



(A4) 



where the factor of 2 accounts for the spin-1/2 degener- 
acy, and the second-order correction. 



5E/A can be obtained by accounting for the density of 
states per area and per spin l/(27r£^) at the unperturbed 
energy -Aq, i.e.. 



5^ 
~A 

If we define the stiffness by 

5E _ 



we deduce that 



27r 



(A6) 



(A7) 



(A8) 



Because the self-consistent Kekulc gap Aq induced by a 
magnetic field and an electron-phonon interaction scales 
linearly in B, the spin-stiffness also scales linearly with 
B. 



-(2) 



IgpAo, (A5) 



are independent of k. The total energy shift per unit area 



Appendix B: Zero modes 



J 



The zero modes of the Dirac Hamiltonian ( ^.^ ), when the Kekulc order parameter A is defective in that it carries 
a vortex of vorticity ±1, can be found analytically and shown to exist for all strengths of the applied magnetic field. 
In this appendix, we choose to work in the symmetric gauge A = ^ {—y,x), which we write as the complex number 
A ~ iBz/2 with z = x -\- iy, and seek the solutions (zero modes) to 





hv^ (-20J - ^ 

A 






A 



A 






A 







* = 0. 



(Bl) 



The spinor 




(B2a) 



is a zero mode supported on sublattice A if 



Au + hvv I 2id^ 



u + Av^O, 



izll 



(B2b) 
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and u and v are normalizable. Switching to polar coordinates and considering an anti- vortex (vorticity 
concreteness, i.e., A = A(p, 9) = Aq e'*" e~*^, the conditions on the components u and v become 



211 



An e-'^^u + ihv^ I d„ + -dg 



21% 



V = 0. 



-1) for 

(B3a) 
(B3b) 



The homogeneous system (B3) of first-order pa rtial differential equations admits ^-independent solutions u(p) and 
v{p). Indeed, solving for u in terms of v in Eq. (B3b) gives 



and, after substitution into Eq. (B3a 



21% 



(B4a) 



(B4b) 



where 



P 



2A^ 



(B4c) 



Equation ( B4h| ) would be identical to the second-order differential equation describing a ID quantum harmonic 
oscillator were it not for the fact that the radial coordinate p is always positive. Equation (B4t) thus admits 
normalizable solutions unavailable to the ID quantum harmonic oscillator, i .e., a ny normalizable solution on the 
half- line that blows up-, when p — > —oo is here allowed. The solutions of Eq. ( B4b ) for p > are parabolic cylinder 
functions D_p{p/ Ib)^ parametrized by the dimensionless ratio p. We conclud e, w hen the Kekule order parameter 
admits an anti- vortex of charg e = —1, that a normalizable zero mode of Eq. (Bl) can be expressed in terms u and 
V satisfying Eqs. (2.21a) and ( ^.21b ). The solution for the case of vorticity ~ +1 is analogous, but with a zero 
mode supported on sublattice B. 
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